Ssid SH Nutrients
Ana Palacio
2020-04-12
General project set-up
# Get all libraries and sources required to run the script
library(plyr)
library(dplyr)
library(reshape2)
library(ggplot2)
library(ggthemes)
library(scales)
library(lme4)
library(lmerTest)
library(emmeans)
library(skimr)
# Default ggplot settings
Fill.colour<-scale_colour_manual(values = c("black", "gray70", "gray35"))
ggthe_bw<-theme(plot.background=element_blank(),
panel.grid.major.y = element_blank(),
panel.grid.major.x = element_blank(),
panel.grid.minor.x = element_blank(),
panel.grid.minor.y = element_blank(),
legend.box.background = element_rect(),
panel.background =element_rect(fill = NA, color = "black")
)+
theme_bw()1. DATA Exploration
1. Select qPCR info and define factors
# Import the data
qPCR.variables <- read.csv("Outputs/Ssid_SH_cell_ratio.csv", header = T)
# Remove blastate samples
qPCR.blastate<-qPCR.variables[((qPCR.variables$Date=="2018-02-05") |
(qPCR.variables$Date=="2018-03-08")), ]
qPCR.variables <- droplevels(qPCR.variables[!rownames(qPCR.variables) %in% rownames(qPCR.blastate), ])
# Variable types
str(qPCR.variables)## 'data.frame': 644 obs. of 17 variables:
## $ Treatment : Factor w/ 3 levels "A","N","N+P": 1 1 1 1 1 1 1 1 1 1 ...
## $ Replicate : Factor w/ 2 levels "R1","R2": 1 1 1 1 1 1 1 1 1 1 ...
## $ Fragment : Factor w/ 162 levels "Ss_20-1","Ss_20-11",..: 27 30 31 48 49 52 55 58 97 100 ...
## $ Genotype : Factor w/ 7 levels "Ss_20","Ss_22",..: 2 2 2 2 3 3 3 3 5 5 ...
## $ Date : Factor w/ 5 levels "2017-11-15","2017-12-14",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ Days : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Time_Point: Factor w/ 5 levels "T13","T19","T22",..: 4 4 4 4 4 4 4 4 4 4 ...
## $ Phase : Factor w/ 3 levels "Baseline","Heat",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ Sample : Factor w/ 644 levels "Ss_20-1_T5","Ss_20-1_T9",..: 109 124 128 188 191 204 217 232 380 393 ...
## $ C.Ssid : num 0 0 0 0 0 ...
## $ D.Ssid : num 0.0266 0.034 0.0239 0.0281 1.0265 ...
## $ TotalSH : num 0.0266 0.034 0.0239 0.0281 1.0265 ...
## $ logC.SH : num NA NA NA NA NA ...
## $ logD.SH : num -1.5756 -1.4689 -1.621 -1.5515 0.0114 ...
## $ logSH : num -1.5756 -1.4689 -1.621 -1.5515 0.0114 ...
## $ D.Prp : num 1 1 1 1 1 ...
## $ Community : Factor w/ 3 levels "C1","C3","D": 3 3 3 3 3 3 3 3 3 3 ... qPCR.variables$Genotype<-factor(qPCR.variables$Genotype,
levels=c("Ss_22","Ss_23","Ss_27", "Ss_28",
"Ss_20", "Ss_24","Ss_30"))
qPCR.variables$DaysF<-as.factor(qPCR.variables$Days)
summary(qPCR.variables)## Treatment Replicate Fragment Genotype Date
## A :220 R1:333 Ss_20-17: 5 Ss_22:92 2017-11-15:140
## N :228 R2:311 Ss_20-19: 5 Ss_23:96 2017-12-14:154
## N+P:196 Ss_20-22: 5 Ss_27:96 2018-02-01:151
## Ss_20-26: 5 Ss_28:82 2018-02-23:103
## Ss_20-29: 5 Ss_20:97 2018-03-06: 96
## Ss_20-32: 5 Ss_24:85
## (Other) :614 Ss_30:96
## Days Time_Point Phase Sample
## Min. : 0.00 T13:151 Baseline :140 Ss_20-1_T5 : 1
## 1st Qu.: 29.00 T19:103 Heat :199 Ss_20-1_T9 : 1
## Median : 78.00 T22: 96 Nutrients:305 Ss_20-11_T13: 1
## Mean : 57.76 T5 :140 Ss_20-11_T19: 1
## 3rd Qu.:100.00 T9 :154 Ss_20-11_T22: 1
## Max. :111.00 Ss_20-11_T9 : 1
## (Other) :638
## C.Ssid D.Ssid TotalSH
## Min. :0.000000 Min. :0.000000 Min. :0.0000034
## 1st Qu.:0.000000 1st Qu.:0.000000 1st Qu.:0.0093391
## Median :0.002373 Median :0.006613 Median :0.0313298
## Mean :0.046419 Mean :0.059660 Mean :0.1060796
## 3rd Qu.:0.027225 3rd Qu.:0.044028 3rd Qu.:0.1286860
## Max. :1.155796 Max. :1.394495 Max. :1.3944953
##
## logC.SH logD.SH logSH D.Prp
## Min. :-5.47264 Min. :-4.7761 Min. :-5.4726 Min. :0.0000
## 1st Qu.:-2.80684 1st Qu.:-2.2410 1st Qu.:-2.0297 1st Qu.:0.0000
## Median :-2.09749 Median :-1.6835 Median :-1.5040 Median :0.6792
## Mean :-2.08904 Mean :-1.7359 Mean :-1.5221 Mean :0.5242
## 3rd Qu.:-1.27217 3rd Qu.:-1.0631 3rd Qu.:-0.8905 3rd Qu.:1.0000
## Max. : 0.06288 Max. : 0.1444 Max. : 0.1444 Max. :1.0000
## NA's :184 NA's :205
## Community DaysF
## C1:132 0 :140
## C3:171 29 :154
## D :341 78 :151
## 100:103
## 111: 96
##
## skim(qPCR.variables)| Name | qPCR.variables |
| Number of rows | 644 |
| Number of columns | 18 |
| _______________________ | |
| Column type frequency: | |
| factor | 10 |
| numeric | 8 |
| ________________________ | |
| Group variables | None |
Variable type: factor
| skim_variable | n_missing | complete_rate | ordered | n_unique | top_counts |
|---|---|---|---|---|---|
| Treatment | 0 | 1 | FALSE | 3 | N: 228, A: 220, N+P: 196 |
| Replicate | 0 | 1 | FALSE | 2 | R1: 333, R2: 311 |
| Fragment | 0 | 1 | FALSE | 162 | Ss_: 5, Ss_: 5, Ss_: 5, Ss_: 5 |
| Genotype | 0 | 1 | FALSE | 7 | Ss_: 97, Ss_: 96, Ss_: 96, Ss_: 96 |
| Date | 0 | 1 | FALSE | 5 | 201: 154, 201: 151, 201: 140, 201: 103 |
| Time_Point | 0 | 1 | FALSE | 5 | T9: 154, T13: 151, T5: 140, T19: 103 |
| Phase | 0 | 1 | FALSE | 3 | Nut: 305, Hea: 199, Bas: 140 |
| Sample | 0 | 1 | FALSE | 644 | Ss_: 1, Ss_: 1, Ss_: 1, Ss_: 1 |
| Community | 0 | 1 | FALSE | 3 | D: 341, C3: 171, C1: 132 |
| DaysF | 0 | 1 | FALSE | 5 | 29: 154, 78: 151, 0: 140, 100: 103 |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| Days | 0 | 1.00 | 57.76 | 41.59 | 0.00 | 29.00 | 78.00 | 100.00 | 111.00 | ▆▆▁▆▇ |
| C.Ssid | 0 | 1.00 | 0.05 | 0.12 | 0.00 | 0.00 | 0.00 | 0.03 | 1.16 | ▇▁▁▁▁ |
| D.Ssid | 0 | 1.00 | 0.06 | 0.14 | 0.00 | 0.00 | 0.01 | 0.04 | 1.39 | ▇▁▁▁▁ |
| TotalSH | 0 | 1.00 | 0.11 | 0.17 | 0.00 | 0.01 | 0.03 | 0.13 | 1.39 | ▇▁▁▁▁ |
| logC.SH | 184 | 0.71 | -2.09 | 1.08 | -5.47 | -2.81 | -2.10 | -1.27 | 0.06 | ▁▃▇▇▅ |
| logD.SH | 205 | 0.68 | -1.74 | 0.94 | -4.78 | -2.24 | -1.68 | -1.06 | 0.14 | ▁▂▆▇▃ |
| logSH | 0 | 1.00 | -1.52 | 0.83 | -5.47 | -2.03 | -1.50 | -0.89 | 0.14 | ▁▁▃▇▅ |
| D.Prp | 0 | 1.00 | 0.52 | 0.47 | 0.00 | 0.00 | 0.68 | 1.00 | 1.00 | ▇▁▁▁▇ |
2. Exploratory graphs
- Histograms
HistoL_SH<-qplot(logSH, data=qPCR.variables, binwidth=0.15)
HistoL_SH + facet_grid(Treatment~Date)
ggplot(qPCR.variables, aes(logSH, fill = Treatment , colour = Treatment)) +
geom_density(alpha = 0.1) + facet_wrap(~Date) + ggthe_bw
- Log SH
logSHGenotype <- ggplot(qPCR.variables, aes(Days, logSH)) +
geom_line(aes(colour=Fragment))+geom_point(aes(shape=factor(Community), colour=factor(Community)))+
# geom_jitter(aes(colour=factor(Replicate))) +
# stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2)+
#stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5) +
# stat_summary(fun.y=mean, geom="line", size =1, alpha=0.5) +
facet_grid(Treatment~Genotype) +
ggthe_bw +theme(legend.position = "none" )
logSHGenotype + ylab("Relative log10 (S:H)") + xlab("Treatment") +
theme(axis.title.y=element_text(size=12), legend.position="none")
logSHTreatment <- ggplot(qPCR.variables, aes(Date, logSH, colour=factor(Treatment))) +
# geom_jitter(aes(colour=factor(Treatment))) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2)+
stat_summary(fun.y=mean, geom="line", size =1) +
theme_bw()
logSHTreatment + facet_wrap(~Genotype)
logSH_Replicate<- ggplot(qPCR.variables, aes (Days, logSH,
colour=factor(Replicate))) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2 )+
stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5) +
stat_summary(fun.y=mean, geom="line") + facet_grid (~Treatment) + ggthe_bw
logSH_Replicate 
logSH_Replicate+ facet_grid(~Community)
logSHTreatment <- ggplot (qPCR.variables, aes(Days, logSH, colour=Treatment)) +
ggtitle("A.") + # geom_point(alpha=0.3) +
ggthe_bw + theme(legend.position="bottom") +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
#stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
# position = position_dodge(width = 1)) +
stat_summary(fun.y=mean, geom="line", position = position_dodge(width = 5)) +
scale_y_continuous(breaks = seq(-5, 0.5, 0.5),
expand = c(0,0),
name=("log 10 (Total S/H cell ratio")) +
scale_x_continuous(name="Days in the experiment",
limits = c(-2,113),
breaks = seq(0, 110, 15),
expand = c(0, 0)) +
annotate("segment", x = 2, xend = 91, y = -3.5, yend = -3.5,
colour = "gray90")+
annotate("segment", x = 79, xend = 90, y = -3.5, yend = -3,
colour = "gray90")+
annotate("segment", x = 91, xend = 110, y = -3, yend = -3,
colour = "gray90")
logSHTreatment
#logSHTreatment + facet_wrap(~Genotype)
logSHTreatment + facet_wrap(~Community)
- Log C
logC_Treatment <- ggplot(qPCR.variables,
aes(Days, logC.SH, colour=factor(Fragment))) +
geom_line()+
#geom_jitter(aes(colour=factor(Replicate))) +
#stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2)+
#stat_summary(fun.y=mean, geom="line", size =1, alpha=0.5) +
ggthe_bw + theme(legend.position = "none" )
logC_Treatment
logC_Treatment + facet_grid(Genotype~Treatment)
logC_Treatment + facet_grid(Community~Treatment)
- Log D
logD_genotype <- ggplot(qPCR.variables,
aes(Days, logD.SH, colour=factor(Fragment))) +
geom_line()+
#geom_jitter(aes(colour=factor(Replicate))) +
#stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2)+
#stat_summary(fun.y=mean, geom="line", size =1, alpha=0.5) +
ggthe_bw + theme(legend.position = "none" )
logD_genotype
logD_genotype + facet_grid(Genotype~Treatment)
logD_Treat <- ggplot(qPCR.variables,
aes(Days, logD.SH, colour=factor(Treatment))) +
#geom_jitter(aes(colour=factor(Replicate))) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2)+
stat_summary(fun.y=mean, geom="line", size =1, alpha=0.5) +
ggthe_bw + theme(legend.position = "none" )
logD_Treat
logD_Treat + facet_grid(Community~Treatment)
*D proportion
D.PTreatment <- ggplot(qPCR.variables, aes(Days, D.Prp, colour=Fragment)) +
ggthe_bw + theme(legend.position="none") +
geom_jitter(aes(shape=factor(Replicate)), alpha=0.5) +
geom_line() + facet_grid (Genotype~Treatment) +
annotate("segment", x = 2, xend = 91, y = 0, yend = 0,
colour = "gray90")+
#annotate("text", x = 45, y = 0.2, label = "Nutrients", size=3)+
annotate("segment", x = 79, xend = 90, y = 0, yend = 0.3,
colour = "gray90")+
annotate("segment", x = 91, xend = 110, y = 0.3, yend = 0.3,
colour = "gray90")+
#annotate("text", x = 99, y = 0.3, label = "Heat", size=3)+
scale_y_continuous(breaks = seq(0, 1, 0.3),
name=("Durusdinium proportion (D/H)/(S/H)")) +
scale_x_continuous(name="Days in the experiment",
limits = c(-2,113),
breaks = seq(0, 110, 45),
expand = c(0, 0))
D.PTreatment
2. Subset data
# Removing weird points and colonies not used
qPCR.variables_2<-subset(qPCR.variables, Date!="2017-12-14")
qPCR.variables_2<-subset(qPCR.variables_2, TotalSH<0.8)
SH.0<-subset(qPCR.variables_2, Days<2)
SH.C<-subset(qPCR.variables_2, Days<79)
SH.C<-subset(SH.C, Days>77)
SH.H<-subset(qPCR.variables_2, Days>80)
logSH_Genotype<- ggplot(qPCR.variables_2, aes (Days, logSH,
colour=factor(Fragment))) +
theme_bw() + geom_line()+ facet_grid (Genotype~Treatment)#+
logSH_Genotype + theme(legend.position = "none" ) + geom_point()
Log S/H graph
logSHTreatment <- ggplot (qPCR.variables_2, aes(Days, logSH, colour=Treatment)) + theme_bw() + ggtitle("A.")+ Fill.colour+ ggthe_bw+
theme(plot.background=element_blank(),
panel.border = element_blank(),
axis.line = element_line(colour = "black"),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position="bottom",
strip.background = element_rect(fill="white")) +
#geom_jitter(aes(colour=factor(Replicate))) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
#stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
# position = position_dodge(width = 5)) +
stat_summary(fun.y=mean, geom="line") +
annotate("segment", x = 2, xend = 91, y = -4.5, yend = -4.5,
colour = "gray90", linetype=1)+
annotate("segment", x = 79, xend = 91, y = -4.5, yend = -4,
colour = "gray90", linetype=1)+
annotate("segment", x = 91, xend = 110, y = -4, yend = -4,
colour = "gray90", linetype=1)+
annotate("text", x = 45, y = -4.4, label = "Nutrients", size=3)+
annotate("text", x = 99, y = -3.9, label = "Heat", size=3) +
scale_y_continuous(breaks = seq(-5, 0.3, 0.5),
expand = c(0,0),
name=("Symbiodiniaceae / S. siderea (S/H) cell ratio")) +
scale_x_continuous(name="Days in the experiment",
limits = c(-2,113),
breaks = seq(0, 110, 15),
expand = c(0, 0))
logSHTreatment 
logSHTreatment + facet_grid(~Community)
#ggsave(file="5.3A_AH.svg", plot=logSHTreatment, width=3.0, height=6)
#ggsave(file="5.3A_AH_NoGenotype.svg", plot=logSHTreatment, width=3.0, height=3)Log C/H graph
CHTreatment_2 <- ggplot (qPCR.variables_2, aes(Days, logC.SH, colour=Treatment)) + theme_bw() + ggtitle("C.")+ Fill.colour+ ggthe_bw+
annotate("segment", x = 2, xend = 91, y = -4.5, yend = -4.5,
colour = "gray90", linetype=1)+
annotate("segment", x = 79, xend = 91, y = -4.5, yend = -4,
colour = "gray90", linetype=1)+
annotate("segment", x = 91, xend = 110, y = -4, yend = -4,
colour = "gray90", linetype=1)+
annotate("text", x = 45, y = -4.4, label = "Nutrients", size=3)+
annotate("text", x = 99, y = -3.9, label = "Heat", size=3) +
#geom_jitter(aes(colour=factor(Replicate))) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
#stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
# position = position_dodge(width = 5)) +
stat_summary(fun.y=mean, geom="line") +
scale_x_continuous(name="Days in the experiment",
limits = c(-2,113),
breaks = seq(0, 110, 15),
expand = c(0, 0))
CHTreatment_2 
#CHTreatment_2 + facet_grid (Genotype~Community)
CHTreatment_2 + facet_grid(~Community)
#ggsave(file="5.3A_AH_2.svg", plot=SHTreatment_2, width=3.0, height=6)
#ggsave(file="5.3C_SH_NoGenotype_2.svg", plot=SHTreatment_2, width=3.0, height=3)Log D/H graph
DHTreatment_2 <- ggplot (qPCR.variables_2, aes(Days, logD.SH, colour=Treatment)) + theme_bw() + ggtitle("C.")+ Fill.colour+ ggthe_bw+
annotate("segment", x = 2, xend = 91, y = -4.5, yend = -4.5,
colour = "gray90", linetype=1)+
annotate("segment", x = 79, xend = 91, y = -4.5, yend = -4,
colour = "gray90", linetype=1)+
annotate("segment", x = 91, xend = 110, y = -4, yend = -4,
colour = "gray90", linetype=1)+
annotate("text", x = 45, y = -4.4, label = "Nutrients", size=3)+
annotate("text", x = 99, y = -3.9, label = "Heat", size=3) +
#geom_jitter(aes(colour=factor(Replicate))) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
#stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
# position = position_dodge(width = 5)) +
stat_summary(fun.y=mean, geom="line") +
scale_x_continuous(name="Days in the experiment",
limits = c(-2,113),
breaks = seq(0, 110, 15),
expand = c(0, 0))
DHTreatment_2 + facet_grid (Genotype~Community)
DHTreatment_2 + facet_grid(~Community)
D Proportion
D.PTreatment <- ggplot(qPCR.variables, aes(Days, D.Prp, colour=Fragment)) +
ggthe_bw + theme(legend.position="none") +
geom_jitter(aes(shape=factor(Treatment)), alpha=0.5) +
geom_line()+
annotate("segment", x = 2, xend = 91, y = 0, yend = 0,
colour = "gray90")+
#annotate("text", x = 45, y = 0.2, label = "Nutrients", size=3)+
annotate("segment", x = 79, xend = 90, y = 0, yend = 0.3,
colour = "gray90")+
annotate("segment", x = 91, xend = 110, y = 0.3, yend = 0.3,
colour = "gray90")+
#annotate("text", x = 99, y = 0.3, label = "Heat", size=3)+
scale_y_continuous(breaks = seq(0, 1, 0.3),
name=("Durusdinium proportion (D/H)/(S/H)")) +
scale_x_continuous(name="Days in the experiment",
limits = c(-2,113),
breaks = seq(0, 110, 45),
expand = c(0, 0))
D.PTreatment + facet_grid (Genotype~Replicate)
D.PTreatment
3. GLM models
Baseline
LME_BaseLine<-lmer(logSH ~ Treatment + Community + (1|Genotype) + (1|Replicate),
data=SH.0)
step(LME_BaseLine)## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 8 -75.016 166.03
## (1 | Genotype) 0 7 -82.009 178.02 13.986 1 0.0001842 ***
## (1 | Replicate) 0 7 -109.286 232.57 68.539 1 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Treatment 1 0.59627 0.29813 2 126.566 2.0641 0.1312
## Community 2 0.79889 0.39945 2 6.535 2.7209 0.1381
##
## Model found:
## logSH ~ (1 | Genotype) + (1 | Replicate) anova(LME_BaseLine) # Treatments is not significant ranova(LME_BaseLine) # Treatments is significant! LME_BaseLine<-lmer(logSH ~ Genotype + (1|Replicate),
data=SH.0)
step(LME_BaseLine)## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 9 -72.767 163.53
## (1 | Replicate) 0 8 -106.895 229.79 68.256 1 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Genotype 0 8.6374 1.4396 6 129 9.7555 7.505e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Genotype + (1 | Replicate) anova(LME_BaseLine) # Treatments is not significant ranova(LME_BaseLine) # Treatments is significant!# Multicomp
Ssid.SH.C.emm<-emmeans(LME_BaseLine, ~Genotype)
#contrast(Ofav.SH.C.emm, "tukey")
Ssid.SH.C_groups<-cld(Ssid.SH.C.emm, by=NULL) # compact-letter display
Ssid.SH.C_groups# Effect plot
plot(emmeans(LME_BaseLine, ~Genotype), comparisons = TRUE) +
coord_flip(xlim = NULL, ylim = NULL, expand = TRUE) + theme_bw()
logSHTreatmentBL <- ggplot (SH.0, aes(Genotype, logSH, colour=Genotype)) +
#ggtitle("A.") +
ggthe_bw + theme(legend.position="bottom") +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
position = position_dodge(width = 1)) +
scale_y_continuous(breaks = seq(-2.5, 0, 0.5),
expand = c(0.1,0.1),
name=("log 10 (Total S/H cell ratio"))
logSHTreatmentBL
logSHTreatmentBL + facet_grid(~Community)
C: Effect of nutrient treatments at control temperature
LME_Ssid.C<-lmer(logSH ~ Treatment + Community + (1|Replicate) + (1|Genotype),
data=SH.C)
step (LME_Ssid.C) # Remove Genotype## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 8 -73.147 162.29
## (1 | Genotype) 1 7 -73.264 160.53 0.235 1 0.6278
## (1 | Replicate) 0 6 -105.884 223.77 65.240 1 6.63e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Treatment 1 0.7873 0.3936 2 144 2.8627 0.06037 .
## Community 0 9.4779 4.7390 2 146 33.6053 9.879e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Community + (1 | Replicate) anova (LME_Ssid.C) ranova (LME_Ssid.C)LME_Ssid.C<-lmer(logSH ~ Treatment + Community + (1|Replicate),
data=SH.C)
step (LME_Ssid.C) # Remove replicate## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 7 -73.264 160.53
## (1 | Replicate) 0 6 -105.884 223.77 65.24 1 6.63e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Treatment 1 0.7873 0.3936 2 144 2.8627 0.06037 .
## Community 0 9.4779 4.7390 2 146 33.6053 9.879e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Community + (1 | Replicate) anova (LME_Ssid.C) ranova (LME_Ssid.C)# Multicomp
Ssid.SH.C.emm<-emmeans(LME_Ssid.C, ~Community)
#contrast(Ssid.SH.C.emm, "tukey")
Ssid.SH.C_groups<-cld(Ssid.SH.C.emm, by=NULL) # compact-letter display
Ssid.SH.C_groups# Effect plot
plot(emmeans(LME_Ssid.C, ~Treatment), comparisons = TRUE) +
coord_flip(xlim = NULL, ylim = NULL, expand = TRUE) + theme_bw()
logSHTreatmentC <- ggplot (SH.C, aes(Treatment, logSH, colour=Treatment, shape=Treatment)) +
ggtitle("A.")+ Fill.colour +
ggthe_bw + theme(legend.position="bottom") +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
position = position_dodge(width = 1)) +
scale_y_continuous(breaks = seq(-2.5, 0, 0.5),
expand = c(0.1,0.1),
name=("log 10 (Total S/H cell ratio"))
logSHTreatmentC
logSHGenotypeC <- ggplot (SH.C, aes(Genotype, logSH, colour=Genotype)) +
ggthe_bw + theme(legend.position="bottom") +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
position = position_dodge(width = 1)) +
scale_y_continuous(breaks = seq(-2.5, 0, 0.5),
expand = c(0.1,0.1),
name=("log 10 (Total S/H cell ratio"))
logSHGenotypeC
H: Effect of pre-exposure to nutrient treatments during heat challenge
LME_Ssid.H<-lmer(logSH ~ Treatment + Community + DaysF+ (1|Replicate) + (1|Genotype),
data=SH.H)
step (LME_Ssid.H) # Remove Replicate## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 9 -150.52 319.04
## (1 | Replicate) 1 8 -150.52 317.04 0.000 1 0.9999
## (1 | Genotype) 0 7 -178.47 370.93 55.896 1 7.641e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Treatment 0 2.8371 1.4186 2 186.93 6.4838 0.001895 **
## Community 0 5.8196 2.9098 2 189.17 13.2996 3.942e-06 ***
## DaysF 0 9.3166 9.3166 1 187.32 42.5827 6.161e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Treatment + Community + DaysF + (1 | Genotype) anova (LME_Ssid.H) ranova (LME_Ssid.H)LME_Ssid.H<-lmer(logSH ~ Treatment + Community + DaysF+ (1|Genotype),
data=SH.H)
step (LME_Ssid.H) # Remove replicate## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 8 -150.52 317.04
## (1 | Genotype) 0 7 -178.47 370.93 55.896 1 7.641e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Treatment 0 2.8371 1.4186 2 186.93 6.4838 0.001895 **
## Community 0 5.8196 2.9098 2 189.17 13.2996 3.942e-06 ***
## DaysF 0 9.3166 9.3166 1 187.32 42.5827 6.161e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Treatment + Community + DaysF + (1 | Genotype) anova (LME_Ssid.H) ranova (LME_Ssid.H)# Multicomp
Ssid.SH.H.emm<-emmeans(LME_Ssid.H, ~Treatment)
#contrast(Ssid.SH.H.emm, "tukey")
Ssid.SH.H_groups<-cld(Ssid.SH.H.emm, by=NULL) # compact-letter display
Ssid.SH.H_groups# Effect plot
plot(emmeans(LME_Ssid.H, ~Treatment), comparisons = TRUE) +
coord_flip(xlim = NULL, ylim = NULL, expand = TRUE) + theme_bw()
logSHTreatmentH <- ggplot (SH.H, aes(Treatment, logSH, colour=Treatment, shape=Treatment)) +
ggtitle("A.")+ Fill.colour +
ggthe_bw + theme(legend.position="bottom") +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar",
position = position_dodge(width = 5), width = 0.2 )+
stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5,
position = position_dodge(width = 1)) +
scale_y_continuous(breaks = seq(-3, 0, 0.5),
expand = c(0.1,0.1),
name=("log 10 (Total S/H cell ratio"))
logSHTreatmentH
logSHTreatmentH + facet_grid(~Community)
logSHTreatmentH + facet_grid(DaysF~Community) #+ geom_jitter()
logSHGenotypeH <- ggplot (SH.H, aes(Community, logSH, colour=Genotype)) +
ggthe_bw + theme(legend.position="bottom") +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2 )+
stat_summary(fun.y=mean, geom="point", size =3, alpha=0.5) +
scale_y_continuous(breaks = seq(-2.5, 0, 0.5),
expand = c(0.1,0.1),
name=("log 10 (Total S/H cell ratio"))
logSHGenotypeH
logSHGenotypeH + facet_grid(~DaysF)
All phases
- Dasys as a factor
LM_Ssid.qpCR<-lmer(logSH ~ Treatment * DaysF * Community + (1|Replicate) + (1|Genotype/Fragment),
data= qPCR.variables_2)
step(LM_Ssid.qpCR)## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df
## <none> 40 -345.92 771.83
## (1 | Replicate) 1 39 -345.92 769.83 0.000 1
## (1 | Fragment:Genotype) 2 38 -345.92 767.83 0.000 1
## (1 | Genotype) 0 37 -376.42 826.84 61.013 1
## Pr(>Chisq)
## <none>
## (1 | Replicate) 1
## (1 | Fragment:Genotype) 1
## (1 | Genotype) 5.669e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value
## Treatment:DaysF:Community 1 2.3178 0.1932 12 443.35 0.9151
## Treatment:Community 2 0.7151 0.1788 4 455.41 0.8489
## Treatment:DaysF 0 5.3896 0.8983 6 459.60 4.2710
## DaysF:Community 0 25.2283 4.2047 6 459.70 19.9919
## Pr(>F)
## Treatment:DaysF:Community 0.531574
## Treatment:Community 0.494734
## Treatment:DaysF 0.000335 ***
## DaysF:Community < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Treatment + DaysF + Community + (1 | Genotype) + Treatment:DaysF +
## DaysF:Community anova(LM_Ssid.qpCR) ranova(LM_Ssid.qpCR)LM_Ssid.qpCR<-lmer(logSH ~ Treatment + DaysF + Community +
(1 | Genotype) +
Treatment:DaysF + DaysF:Community,
data= qPCR.variables_2)
step(LM_Ssid.qpCR)## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 22 -344.90 733.79
## (1 | Genotype) 0 21 -376.15 794.30 62.504 1 2.659e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Treatment:DaysF 0 5.3896 0.8983 6 459.6 4.271 0.000335
## DaysF:Community 0 25.2283 4.2047 6 459.7 19.992 < 2.2e-16
##
## Treatment:DaysF ***
## DaysF:Community ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Treatment + DaysF + Community + (1 | Genotype) + Treatment:DaysF +
## DaysF:Community anova(LM_Ssid.qpCR) ranova(LM_Ssid.qpCR) summary(LM_Ssid.qpCR)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## logSH ~ Treatment + DaysF + Community + (1 | Genotype) + Treatment:DaysF +
## DaysF:Community
## Data: qPCR.variables_2
##
## REML criterion at convergence: 689.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3945 -0.6829 -0.0445 0.7223 3.3266
##
## Random effects:
## Groups Name Variance Std.Dev.
## Genotype (Intercept) 0.4366 0.6608
## Residual 0.2103 0.4586
## Number of obs: 486, groups: Genotype, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.91537 0.28599 8.61510 -6.697 0.000109 ***
## TreatmentN -0.07283 0.09326 459.24915 -0.781 0.435249
## TreatmentN+P -0.23402 0.09857 459.28241 -2.374 0.018003 *
## DaysF78 0.18484 0.14532 459.31173 1.272 0.204047
## DaysF100 -0.87819 0.15269 459.54363 -5.752 1.62e-08 ***
## DaysF111 -1.16397 0.15636 459.47375 -7.444 4.85e-13 ***
## CommunityC3 1.46403 0.22145 395.52840 6.611 1.24e-10 ***
## CommunityD 0.16239 0.14456 465.99482 1.123 0.261877
## TreatmentN:DaysF78 0.09794 0.13091 459.24978 0.748 0.454730
## TreatmentN+P:DaysF78 0.08522 0.13522 459.27412 0.630 0.528863
## TreatmentN:DaysF100 0.48788 0.14188 459.32643 3.439 0.000638 ***
## TreatmentN+P:DaysF100 0.43816 0.15142 459.30234 2.894 0.003990 **
## TreatmentN:DaysF111 0.09205 0.14651 460.63193 0.628 0.530117
## TreatmentN+P:DaysF111 0.48196 0.15444 459.53356 3.121 0.001917 **
## DaysF78:CommunityC3 0.35575 0.16130 459.31491 2.205 0.027917 *
## DaysF100:CommunityC3 -0.12994 0.17289 459.43506 -0.752 0.452701
## DaysF111:CommunityC3 -1.18260 0.19715 459.34228 -5.999 4.04e-09 ***
## DaysF78:CommunityD 0.10776 0.14580 459.32787 0.739 0.460215
## DaysF100:CommunityD 0.02428 0.15623 459.39217 0.155 0.876552
## DaysF111:CommunityD 0.31536 0.15824 460.63181 1.993 0.046859 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Emmeans
Ssid.YII.emm<-emmeans(LM_Ssid.qpCR, ~Treatment*DaysF*Community)
#contrast(Ssid.YII.emm, "tukey")
# Pairwise comparisons
#pairs(Ssid.YII.emm) # same than contrast(Ssid.YII.emm, "tukey")
plot(emmeans(LM_Ssid.qpCR, ~Treatment|DaysF), comparisons = TRUE) +
theme_bw() # Tukey comparission (do not trust CI to compare EMMs)
plot(emmeans(LM_Ssid.qpCR, ~Treatment|DaysF|Community), comparisons = TRUE) +
coord_flip(xlim = NULL, ylim = NULL, expand = TRUE) + facet_grid(Community~DaysF) +
theme_bw()
#CLD
SH_Multicomp<-cld(Ssid.YII.emm, by=NULL) # compact-letter display
SH_Multicomp<-SH_Multicomp[order(SH_Multicomp$DaysF,
SH_Multicomp$Community,
SH_Multicomp$Treatment), ]
write.csv(SH_Multicomp, "Ssid_SH_Multicomp.csv") # 2. Predict values:
pred_Ssid1 <- predict(LM_Ssid.qpCR, re.form = NA)
#3. Bootstrap CI:
Ss.boot1 <- bootMer(LM_Ssid.qpCR, predict, nsim = 1000, re.form = NULL) # include random effects, reduce CI lot!
std.err <- apply(Ss.boot1$t, 2, sd)
CI.lo_1 <- pred_Ssid1 - std.err*1.96
CI.hi_1 <- pred_Ssid1 + std.err*1.96
#Plot
Model_Ss_1b_plot<- ggplot(
qPCR.variables_2, aes(x = Days, y = logSH, colour = Treatment)) +
geom_line(aes(y = pred_Ssid1),size=2) +
#geom_point(aes(fill=factor(Treatment)),
# shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.3) +
geom_ribbon(aes(ymin = CI.lo_1, ymax = CI.hi_1),
size=2, alpha = 0.1, linetype = 0) +
#scale_color_manual(values=my_colours) +
#scale_fill_manual(values=my_colours) +
#scale_y_continuous(name=expression(~italic("Fv / Fm")),
# limits = c(0.5, 0.61),
# breaks = seq(0.4, 0.6, by=0.1), expand = c(0,0))+
#scale_x_continuous("Days in the experiment", limits = c(0, 78),
# breaks = seq(0, 76, by=7), expand = c(0,0))+
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 1,
position = position_dodge(1) )+
stat_summary(fun.y=mean, geom="line", position = position_dodge(1),
linetype=1, alpha=1) +
# stat_summary(fun.y=mean, geom="point", size =1,
# position=position_dodge(width=1), alpha=0.8) +
ggthe_bw
Model_Ss_1b_plot + facet_wrap(~Community)
- Days as continous
LM_Ssid.qpCR2<-lmer(logSH ~ Treatment * Days * Community + (1|Replicate) + (1|Genotype/Fragment),
data= qPCR.variables_2)
step(LM_Ssid.qpCR2)## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df
## <none> 22 -567.31 1178.6
## (1 | Replicate) 1 21 -567.31 1176.6 0.000 1
## (1 | Fragment:Genotype) 2 20 -567.31 1174.6 0.000 1
## (1 | Genotype) 0 19 -588.46 1214.9 42.303 1
## Pr(>Chisq)
## <none>
## (1 | Replicate) 1
## (1 | Fragment:Genotype) 1
## (1 | Genotype) 7.817e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value
## Treatment:Days:Community 1 1.9970 0.4992 4 461.13 1.0674
## Treatment:Community 2 3.0175 0.7544 4 465.28 1.6124
## Days:Community 3 2.0249 1.0125 2 469.82 2.1535
## Community 0 23.5575 11.7787 2 371.44 24.9107
## Treatment:Days 0 3.1263 1.5631 2 471.54 3.3058
## Pr(>F)
## Treatment:Days:Community 0.37206
## Treatment:Community 0.16997
## Days:Community 0.11722
## Community 7.047e-11 ***
## Treatment:Days 0.03752 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Treatment + Days + Community + (1 | Genotype) + Treatment:Days anova(LM_Ssid.qpCR2) ranova(LM_Ssid.qpCR2)LM_Ssid.qpCR2<-lmer(logSH ~Treatment + Days + Community + (1|Genotype) + Treatment:Days,
data= qPCR.variables_2)
step(LM_Ssid.qpCR2)## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 10 -542.25 1104.5
## (1 | Genotype) 0 9 -563.55 1145.1 42.618 1 6.656e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Community 0 23.5575 11.7787 2 371.44 24.9107 7.047e-11
## Treatment:Days 0 3.1263 1.5631 2 471.54 3.3058 0.03752
##
## Community ***
## Treatment:Days *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## logSH ~ Treatment + Days + Community + (1 | Genotype) + Treatment:Days anova(LM_Ssid.qpCR2) ranova(LM_Ssid.qpCR2) summary(LM_Ssid.qpCR2)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## logSH ~ Treatment + Days + Community + (1 | Genotype) + Treatment:Days
## Data: qPCR.variables_2
##
## REML criterion at convergence: 1084.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.5339 -0.6119 -0.1030 0.6531 2.3803
##
## Random effects:
## Groups Name Variance Std.Dev.
## Genotype (Intercept) 0.6894 0.8303
## Residual 0.4728 0.6876
## Number of obs: 486, groups: Genotype, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.889571 0.359416 8.086576 -5.257 0.000741 ***
## TreatmentN -0.097590 0.136063 471.087555 -0.717 0.473579
## TreatmentN+P -0.267215 0.143404 471.032500 -1.863 0.063032 .
## Days -0.008154 0.001207 471.192021 -6.755 4.21e-11 ***
## CommunityC3 1.865263 0.278743 295.864217 6.692 1.10e-10 ***
## CommunityD 0.305446 0.160114 461.670216 1.908 0.057053 .
## TreatmentN:Days 0.003255 0.001699 471.935734 1.916 0.055934 .
## TreatmentN+P:Days 0.004351 0.001797 471.117838 2.421 0.015863 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TrtmnN TrtN+P Days CmmnC3 CmmntD TrtN:D
## TreatmentN -0.183
## TreatmntN+P -0.183 0.480
## Days -0.246 0.591 0.564
## CommunityC3 -0.356 -0.027 0.001 0.051
## CommunityD -0.363 -0.009 0.006 0.044 0.563
## TrtmntN:Dys 0.145 -0.834 -0.400 -0.705 0.070 -0.007
## TrtmntN+P:D 0.137 -0.399 -0.838 -0.669 0.028 0.033 0.477 # 2. Predict values:
pred_Ssid1 <- predict(LM_Ssid.qpCR2,re.form = NA)
#3. Bootstrap CI:
Ss.boot1 <- bootMer(LM_Ssid.qpCR2, predict, nsim = 1000, re.form = NULL) # include random effects, reduce CI lot!
std.err <- apply(Ss.boot1$t, 2, sd)
CI.lo_1 <- pred_Ssid1 - std.err*1.96
CI.hi_1 <- pred_Ssid1 + std.err*1.96
#Plot
Model_Ss_1b_plot<- ggplot(
qPCR.variables_2, aes(x = Days, y = logSH, colour = Treatment)) +
geom_line(aes(y = pred_Ssid1),size=2) +
#geom_point(aes(fill=factor(Treatment)),
# shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.3) +
geom_ribbon(aes(ymin = CI.lo_1, ymax = CI.hi_1),
size=2, alpha = 0.1, linetype = 0) +
#scale_color_manual(values=my_colours) +
#scale_fill_manual(values=my_colours) +
#scale_y_continuous(name=expression(~italic("Fv / Fm")),
# limits = c(0.5, 0.61),
# breaks = seq(0.4, 0.6, by=0.1), expand = c(0,0))+
#scale_x_continuous("Days in the experiment", limits = c(0, 78),
# breaks = seq(0, 76, by=7), expand = c(0,0))+
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 1,
position = position_dodge(1) )+
stat_summary(fun.y=mean, geom="line", position = position_dodge(1),
linetype=1, alpha=1) +
# stat_summary(fun.y=mean, geom="point", size =1,
# position=position_dodge(width=1), alpha=0.8) +
ggthe_bw
Model_Ss_1b_plot + facet_wrap(Genotype~Community)